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Abstract 

Electrical capacitance tomography (ECT) is used to image cross-sections of 
industrial processes containing dielectric material. This technique has been 
under development for more than a decade. The task of image 
reconstruction for ECT is to determine tlie permittivity distribution and 
hence material distribution over the cross-section from capacitance 
measurements. There are three principal difSculdes with unage 
reconstruction for ECT: ( 1 ) the relationship between the permittivity 
distribution and capacitance is non-linear and the electric field is distorted 
by the material present, the so-called 'soft-field* effect; (2) the number of 
independent measurements is limited, leading to an under-determined 
problem and (3) the inverse problem is ill posed and ill conditioned, making 
the solution sensitive to measurement errors and noise. Regularization 
methods are needed to treat this ill-posedness. This paper reviews existing 
image reconstruction algorithms for ECT, including linear back-projection, 
singular value decomposition, Tikhonov regularization, Newton-Raphson, 
iterative Tikhonov, the steepest descent method, Landweber iteration, the 
conjugate, gradient method, algebraic reconstruction techniques, 
simultaneous iterative reconstruction techniques and model-based 
reconstruction. Some of these algorithms are examined by simulation and 
experiment for typical permittivity distributions. Future developments in 
image reconstruction for ECT arc discussed. 

Keywords: electrical capacitance tomography, image reconstruction, iterative 
algoridim, inverse problem 
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List of symbols 

A area of electrode plate 

C. A C capacitance and change in capacitance 

C. AC vectors of capacitance and change in capacitance 

d distance between two electrode plates 

a measurement error vector 

^ Atntior to whom any correspondence shoald be addrcsced. 



/ function 

F transform function from perminivity 

distribution tu capacitance 

g normalized permittivity vector 

g reconstructed permittivity vector 

|. i mean values of g and g respectively 

g estimated solution 

g solution of ^ in null space 
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J identity matrix 

J Jacobian matrix 

k iteration step number 

M number of measurements 

n number of electrodes 

S number of pixels 

p rank of 5 

. p projection operator 

Q electric charge 

s sensitivity of capacitance transducer 

kih row of S 

S normalized sensitivity matrix 

ux Identity normalized capacitance vector 

Um mth column of C;* 

V orthogonal and diagonal matrices for singular 

value decomposition 

v„ rtth column of V 

V potential difference between two elecJtrodes 

wi weighting factor 

X variable vector 

or, p, Y relaxation factors 

S( ith singular value of S 

e(x, y) permittivity distribution 

Bq permittivity of vacuum space 

Br relative permittivity of material 

(Ae)^ second order term in a differential equation 

0{{Asp) {Asf and higher order terms 

<Kx»y) potential distribution 

r electrode surface 

X normalized capaci tance vector 

\t i^th component of A 

Ae noise-contaminated Ae 

A« calculated vector of normalized capacitance 

fi Tlkbonov regularization parameter 

BtL.r model parameters 

I function 



1. Introductioii 

Electrical capacitance tomography (ECT) is based on mea- 
suring capacitances of a mulu-electrode sensor surrounding 
an industrial vessel or pipe containing two materials of differ- 
ent permittivities* These measuremenbs are used to reconstruct 
the permittivity distribution and hence the material distribution 
over the cross section, using a suitable algorithm. Invesliga- 
nons into ECT date back to die 1970s. The first reported ECT 
system was developed by the US Department of Energy in 
Morgantown, to image the gas/solids distribution in fluidized 
beds (Fasching and Smith 1991. Halow et al 1993, Fasching 
et al 1994). The first real-time ECT s>'stem was developed at 
UMtST in die late 1980s to image gas/oil pipelines (Huang 
et al 1988. 1992. Xie ct al 1992). Subsequently, the UMIST 
ECT system has been used successfully in a number of research 
investigations for imaging other two-phase processes, such as 

• gas/solid distribution in pneumatic conveyors and 
fluidized beds (WUliams and Xie 1993, Wang et al 1995, 
Dyakowski etal 1 997^ 2000, Williams ei al 2000, Liu ei al 
2001. iaworski and I^akowski 2001); 

• combustion flames in engine cylinders (He et al 1994, 
WaierfaU 1996); 
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• a water hammer (Yang et al 1 996). 

In the past decade, research into ECT has been carried 
out in several countries (Beck and Williams 1996. Beck 
et al 1997, Hanuner and Johansen 1997. Beck et al 1998, 
Rcineckc and Mcwes 1996, Reinecke et al 1998, York 2001). 
The institutions involved in development of hardware and/or 
software of ECT systems are given In table 1 . 

A typical ECT system as shown in figure 1 comprises diree 
main units: 

(1) a multi-electrode sensor, 

(2) the sensing electronics and 

(3) a computer for hardware control and data processing, 
induding image reconstruction. 

A set of electrodes, usually eight or 12, is mounted around a 
pipe or vessel and the capacitance values between all single- 
electrode combinations are measured. The sensing electronics 
provides excitation signals and converts the capaciumces into 
voltage signals, which are conditioned and then digitized for 
data acquisition. The computer controls the system hardware 
and implements image reconstruction to show the permittivity 
distribution. 

There are two major computational problems in ECT: 
die forward problem and the inverse problem. The forward 
problem is to determine uiter-electrode capacitances from the 
permittivity distribution, e.g. by solving the partial differential 
equations, governing the sensing domain. The inverse problem 
is to determine the permittivity distribution from capacitance 
measurements. The result is usually presented as a visual 
image, and hence this process is called image reconstruction. 
(Current image reconstruction methods, particularly iterative 
methods, also involve solving the forward problem. 

The relationship between capacitance and permittivity 
distribution is governed by the following equation: 

C-| = -:^ jJe{x,y)S;<f>(x.y)6T (1) 
r 

where e{x,y) is die permittivity distribution in the sensing 
field. V is the potential difference between two electrodes 
fonning the capacitance, ^(;c. y) is the potential distribution 
and r is the eIecu*ode surface. 

In some cases, equation (1) can be simplified. For 
example, for an ideal parallel-plate sensor with a homogeneous 
permittivity distribution, it b«:omes 

C:=zeQer- (2) 

where eo is the permittivity of vacuum. Br is the relative 
permiuivity of the material inside the sensor, A is the area 
of the plates and d is die distance between die two plates. 

Equation (2) implies such a simple relationship that 
capacitance Ls proportional to permittivity. However, die 
geometry of an ECT sensor is much more compKcated than 
that of a parallcNplate sensor and the permittivity distribution 
is generally not uniform. In this case, equation (1) cannot 
be simplified. In equation (1), die potential distribution 
^(.v. y) also depends on the permiuivity distribution e(jr. y). 
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Figure 1. Typical ECT system with an eight-cIectrodc sensor, with seasing field boundaries shown. 
Table 1. Institutions involved in development of ECT systems. 



Country Insutution 



Hardware 



Algorithm* 



Applicntion 



Reference 



China Chinese Academy 

of Sciences 
Northeast University 



Tsinghua 
University 
Zhejiang University 
Gcrmaiiy University 
of Haimovcr 



AC-based circuit 
AC-based circuit 



Self-balancing 
AC-based circuit 
Charge transfer 

Capacitance meter/ 
impedance analyser 



LBF, Landweber Fluidized bed 
iteration 

LBP, linear 
regression and 
regularization 

Regularization, 
SIRT 

Oil pipeline 

LBP.ARXSTRT, Trickle bed 
electrical field lines 



Multiphase flow 



Holland Delft University of 

Technology 
Mexico Institute Mexico 

Petroleum 

Norway University of Bergen Self-balancing 

AC'based circuit 



Active difTemntlator 



AC-based circuit LBP 



South Africa Cape Technilcon 
UK Manchester 
Metropolitan 
University 

UMIST 



Self-balancing 
AC-based drcuit 



Charge/dischaxse» 
AC-based circuit, 
LCRmextT 



Model-based 
iteradon 



LBP 



LBP, Landweber 

iteration, 

regularization 



Fluidized bed 

Oil pipeline 

Oil pipeline, 
oil separator 



Gas/solids flow 
metering 

Oil pipeline, 
pneumatic conveyor, 
fluidized bed, 
combustion flame, 
water hammer, 
wet gas 



Liu e/ fly (1999, 2001) 
YanetaH200i) 



Vcng eta( (2000) 
Su era/ (2000) 

Hkhonov and Arsenin (1977) 

Loser e/ a/ (2001) 
Reineclcc and Mewes 
(1994, 1996. 1997) 
Reinedie€fa/(1998) 
Kuhn and van Haideren (1997) 

Gamio (2002) 

Hammer and Jobansen (1997) 
Isalcsen and Noidtvedt (1993) 
Isuksen era/ (1994) 
Isaksen (1996) 

Scott and Gutsche (1999) 

Deloughryefd/ (2001) 



USA 



US Depanmcni 
of Energy, 
Morgantown 



High voltage 
AC bridge 



Fluidized bed 



Beck e/fl/ (1997) 

Dyakowski ei al (1997. 2000) 

Hce/fl/(1994) 

Huang « 0/(1988. 1992) 

Javvonki and £)yakows1ci (2001) 

Uonheart (2001) 

Wang eiff/ (1995) 

Waterfall e; a/ (1996) 

Williams and Xie (1993) 

Williams e/ a/ (2000) 

Xieeia/ (1992. 1994) 

Yang e/ a/ (1995. 1996, 1999) 

Yang (1996) 

Yang and Yorij (1999) 

YorkafidDuggan(1997) 

Yoric(2001) 

Fosching and Smith (1991) 
Fascbing et al (1994) 



* T^e abbreviations are in the main text. 

Therefore, the capacitance between electrode combinations 
can be considered as a functional of permittivity distribudon. 



The change in capacitance in response to a perturbation of the 
permittivi^ distribution is given by 



(3) 



AC=^(A£) + 0((A€)^) 



(4) 
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where g is the sensitivity of the capacitance versus 
permittivity distribution and 0{{^ef) represeats (Ae)^ and 
higher Older tenm. 

Tt) ECT applications, is usually s/nall By neglecting 
0((Ae)^), equation (4) can be simplified to the linear form. 

AC = .vAff (5) 

where 5 =i ^ is the sensitivity of the capacitance transducer 
to changes in permittivity. 

Equation (5) has to be discretized for implementation, in 
order to visualize the permittivity distribution, the sensing area 
is divided into elements or pixels, typically of the order of 
1000. For example, a 32 x 32 grid generates 1024 pixels in 
the image area of a square sensor and about 800 pixels for a 
circular sensor. 

.Using an eighi-eJectrode system as an example, 
capacitance mefisuremenls are obtained in the following steps. 
A voltage sigoal is fini applied lo electrode 1, and then (he 
electric charges on electrodes 2-^ respectively arc measnred, to 
obtain the capacitances between electrode I and the other seven 
electrodes. Next electrodes 2-7 are energized in sequence, so 
that capacitances of all single-electrode pairs are measured, 
giving a total of 28 independent electrode combinations. There 
are = n(n — l)/2 independent capacitance measurements 
for an n-electrode sensor, giving M equations in the form of 
equation (5). The linearized and discrete form of the forward 
problem can now be expressed as 

AC= J Ae (6) 

where J is a Jacobian matrix, i.e. the seasitivity distribution 
matrix, giving a sensitivity map for each electrode pair. 

By this ureatmcDt, the non-linear forward problem has 
been simplified to a linear approximation. This satisfies most 
applications with a small permittivity contrast or pcnurbation. 
Commonly, equation (6) is written in a normalized form (Xie 
etal 1992), 

X=S9 (7) 

where A is the normalized capacitance vector, S is the 
Jacobian matrix of normalized capacitance with respect to 
normalized permittivity, i.e. the transducer sensitivity, and is 
the normalized pennittivity vector, i.e. Ihc grey level of pixels 
for visualization. 

The mauix S contains M sensitivity distributions, but 
with a symmetrical electrode system there will be only n/2 
unique distributions and other disoibutions can be obtained by 
rotating and mirroring. Figure 2 shows two typical sensitivity 
distributions for an adjacent electrode pair and an opposing 
electrode pair respectively, which have been obiained using 
a finite element software package^PC-Opera from Vector 
Fields Ltd (1992), based on a sensitivity theorem (Xie et a! 
J992). 

The task of image reconsUruction for ECT is to determine 
the permittivity distribution e(x,y) from the measured 
capacitance vector C. In the discrete form, it is to find tiie 
unknown g from the known A using equation (7), while S is 
treated as a constant matrix for simplicity. Note that S will 
change with permiltivity distribution. 



(a) 




Figure 2. Sensitivity distributions for an eight-electrode system. 

There are two major difficulties assucimed with 
equation (7). Firstiy, it is under-determined because the 
number of unknown variables N (I.e. the number of pixels) 
is usually much larger than the number of cquatiims M 
(i.e. the number of capacitance measurements). Therefore, the 
solution is not unique. Secondly, equation (1) is an integral 
equation, which is ill posed, and its corresponding discrete 
form, equation (7), is ill conditioned. This means that the 
solution of equation (7) is .sensitive to small perturbations of 
A Cfikhonov and Arsenin 1977, Engl et al 1996). In addition, 
the matrix S is not truly constant, but varies witii die acmal 
permittivity distribution. 

Image reconstruction algorithms for ECT were reviewed 
some years ago (Isalcsen 1996). Since then a number of new 
or different algorithms have been developed to address tiic 
ill posed and ill conditioned problems. In general, -tiiey can 
be categorized into two groups, non-iterative (or single step) 
algoritiims and iterative algorithms. 

2. Non-iterative algoritfams 
2.7. Linear back-projection (LBP) 

If the inverse of S existed, equation (7) would be solved diicctiy 
by 

p = 5-*A. (8) 

However, the inverse of 5 does not exist and so other metiiods 
for solution must be used. The iirst attempt to reconsoruct 
images for ECV was made using linear back-projection (Xie 
ei al 1992). If S is considered lo be a linear mapping 
from the perminivity vector space lo the capacitance vector 
space, can be considered as a related mapping from die 
capacitance vector space to the penrnttixqty vector .space, 
giving an appro-ximaied solution. 

g^S^X. (9) 

A normalized or weighted form of equation (9) can be 
expressed as 
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9 = 



= [1. 1 1] 



(10) 



where is an identity vector. The division of the two 
vectors 5^ A and is defined as one numerator component 
being divided by the corresponding denominator component 
Note that mathematically g is only an approximate solution of 
equation (7). 

If the sensoris filled first with a lower permittivity material 
and Aen with a higher permittivit}' material, the Dormalized 
change In capacitance must be ux and the corresponding 
change in pennittivity is 5^ttx. If the normalized measured 
capacitance is A, the corresponding change in permittivity 
is given by equation (10). Although mathematically not 
accurate, the LBP algorithm is still widely used for on-line 
Im^e reconstrocdon because of its simplicity. However, it 
produces poor-quality images and can only provide qualitative 
information. To obtain improved images or quantitative 
information, other image reconstruction algorithms have to be 
used 

2.2. Direct method based on singular value decomposition 

When measurement errors are considered, equation (7) 
becomes 



Sg = A + e 



(11) 



where e is the capacitance measuiement error vector. 

A common method to obtain the solution of equation (11) 
is to minimize ^||^|^ - A||^, giving 



(12) 



If the inverse of S^S existed, the solution of equation (7) in 
the least squares sense would be 



^ = (S^5)-*5^A. 



(13) 



However, the inverse of 5^5 does not exist in most cases, 

because the number of pixels is usually much larger than 
the number of measurements. Singular value decomposition 
(SVD) provides a means to obtain the solution of equation (7) 
directly. If its rank is p, S can be decomposed as 



where 



S^UTAT'^ (14) 

U = (^1,1^3 um] 

y^[v..v^ vh] (15) 

2:=diagf3|.52,.... 

and U is an M X A/ orthogonal matrix, V is an iV x ^ 
orthogonal matrix (note that any two columns ofU otV are 
orthogonal), S is an ilf x matrix with all components zero 
except for the diagonal components and 52, • . • . ^p-u^p 
are /? singular values of 5 (note that 5j ^ ^2 > ^ Sp^t ^ 
Sp > 0). 

The solution of equation (7) can now be written as 



where Vr~*C/^ is the pseudo-inverse of S, and is given 

The solution of equation (7) or (1 1) is an x ^ diagonal 
matrix. However, it is not unique, because the inverse problem 
is ill posed and under-determined* Suppose that g Is in the null 
space of 5, then 

55 = 0 (18) 

where g can be expressed as a linear combinntion of column 
vector Vp+\, ,Vf^, i.e. the column vector of orthogonal 
matrix V. 

Jfg is the solution of equation (7) or (1 1), g +p must also 
be a solution. The solution according to equations (14)-<16) 
is the minimum norm solution. However, it is only acceptable 
ma±emaUcally. 

Truncated singular value decomposition (TSVD) may 
improve tiiis situation, eKpeciaUy when an obvious gap exists 
between the smgular values. Equation (17) may be modified 
as 

W2 Wp^i 



—I 



Wp 



(19) 



(20) 



where Wi can be considered as a filter and is a rcgularization 
parameter, which must be positive. 

The solution according to equations ( 1 6), ( 1 9) and (20) is 
actuafly die same as Iikhonov rcgularization (to be discussed 
in section 2,3). A simple way of choosing Wf is to let = 1, 
when Sf ^ fi, and wt^Q when 6j < fjL, This is the so-called 
TSVD. Compared to the conventional SVD, TSVD may be 
considered as having a filter, and hence it is less sensitive to 
high frequency noiius in the measurements. 

2,3, Hkhonov regularizfltion 

Regularization tools have been developed for many years to 
solve ill pased inverse problems (Tikhonov and Arscnin 1977, 
Engl et al 1996, Bertero and Boccacci 1998, Hansen 1998). 
They are used to determine a set of solutions using prior 
constraint information and then one solution horn this set 
is cho.sen. Hkhonov method is one of the most commonly 
used universal regularization tools for solving ill posed inverse 
problems andhas been applied toECT for imagereconstruction 
(Peng et al 2000. Lionheart 2001). Based on the standard 
Tikhonov regularization procedure, tiie solution of equation (7) 
or (1 1) is expressed as 



(21) 



whero / is an identity matrix, and /t is the regularization 
parameter, which must be positive. 

Equation (21) may be considered as a modified version of 
equation (13). To make its inveree exist, matrix S^S must be . 
modified into 5^5 + iil. 

A more genend form of Tikhonov regularization is to find 
g, which minimizes the following function: 



9 = VX->C/^A 



(16) 



i(05p-Ay^ + M0X'(p-5)f) 



(22) 
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where ^ is the estimated solution according lo sonie prior 
information and fiiLig — g)i^ is used as a coostraint. 

Minimizing this frnicdon means that the discrepancy 
between the measured capacitance and the esdmated 
capacitance is muiimized while the estimated solution is kept 
reasonably dose to the true solution. Actually, this principle 
underlies most mgularizatlon methods. In pmcdce, it is 
difficuh to obtain g when prior information is not available. In 
genera), choosing g as zeio and li as an identity' matrix gives 
die standard Tikhonov regularizaUon, which is equivalent to 
solving the following equation: 

{S^5 + MDff = 5^A. (23) 

It has been proved that the inverse of 5^5 + ;z/ always exists, 
so equation (21) can be obtained from equation (23). 

The quality of Tikhonov regularizalion strongly depends 
on the regularization parameter ^ It is crucial to choose 
an optimal regularizaUon parameter jx, so thit a solution as 
close to die true solution as possible can be obtained. Some 
mathematical methods, e.g. the discrepancy principle method 
(Engl and Gfrerer 1988), need prior information on the noise 
in the measurements. Other methods, e.g. the generalized 
cross-validation method (Golub et al 1979) and the L-curve 
mediod (Hansen 1992, Hansen and O'Leary 1993), need less 
prior information but the related calculations are laborious. 
In general, a small value of fi gives a good approximation to 
the original problem but the influence of errors may make die 
solution physically unacceptable. Conversely, a large value of 
/z suppresses the data errors but increases die approximation 
error. At present, in almost all cases, fi is chosen empirically. 

2.4, Multiple linear regression and regularization 

A multiple linear regression and regularization algorithm for 
ECT is reported by Yan et al (2001). In principle, it is similar 
to Hkbonov regularization. The main difference is that the 
sensitivity matrix is replaced by a system matrix obtained 
tiirough multi*Iinear regression. 

3« Iterative algorithms 

Because of (he non-h*near relationship between the permittivity 

distribution and the capacitance, it is almost impossible 
to find an accurate solution by any single-step (i.e. non- 
iterative) algorithm with a simplified linear model. To 
improve die image quality, die inverse problem has to be 
solved iteratively. Iterative algorithms developed and applied 
to EOT are based on calculating capadtance values from 
the permittivity distribution of the current image, and then 
producing a new image using the discrepancy between the 
measured capacitance and the calculated capacitance. This 
process is repeated until a satisfactorily low discr^ancy is 
achieved. 

Calculating capacitance values (i.e. solving die forward 
pniblem) is carried out using finite element techniques or by 
superimposing sensitivity maps according to the estimated 
image. This is a time-consuming process and consequentiy 
iterative algorithms are slow compared widi die non-iterative 
aIgorithm.H. Therefore, iterative aigoridims can only be used 
for off-line data processing at present. Dtfferem iterative 



algorithms have been developed, according to the methods 
used to modify die image. 

3 J. Newton-Raphson and tierutive Tikhonov methods 

Newton-Raphson method is an iterative approach, which was 
initially devised to find die root of a non-linear function. For 
image reconstruction widi ECT, suppose J* is die transform 
function from permittivity distribution g to capacitance A; the 
relationship between diem without measurement error can be 
described as 

\^F(g). (24) 

The solution of equation (24) can be deduced by minimizing 
the square error between the measured capacitance and the 
capadtance estimated from die image, i.e. F(g): 

e = l[ii'(y)-AriP(p)-A]. (25) 

Differentiating the above error witii respect lo permittivity 
gives 

e^[F'(g)nF(9)-'X] (26) 

where [F'(g)]ij = |^ is the so-called Jacobian matrix. 

To Rnd die solution of equation (24), let equation (26) 
be zero. Expanding equation (26) in Taylor series and 
neglecting die non-linear or higher order terms give a solution 
of equation (24) as 

9M =9k - [{F'{gx)V{F'{g{))r\{r{g,)f{F{3^ - A)I 

/ ^^^^ 

where F'igk) and F(gii) are the Jacobian and die capacitance 
which are associated with die current pcrmlnivity distribution 
9k' 

Note that in principle F'{3k) £Lnd F{gk) should be updated 
at each iteration step. This is the so-called Newion-Raphson 
method. A simplified form of equation (27) is 

5k! = & - (SfStr^SfiXi ^ A) (28) 

where Sjt is the sensitivity matrix (i.e. Jacobian) and A^ is the 
capadtance. 

However, a problem with equation (28) is that die inverse 
of Sj' Sk usually does not exist because J^S^ is not a full- 
rank matrix. A modified Newton-Raphson iierative algorithm 
can be implemented by adding a term y/ to 5^Sjt. making its 
inverse exist 

PA+i =gk- (5f Si + yl)-^ Sf(\t - A) (29) 

where / is an identity matrix, and y is a positive scalar. 

In optimization theory, dils mediod is also called the 
Levcnburg-Mnnquardt method when y changes with each 
iteration. 

Updating die sensitivity matrix Sk and capacitance A^^ at 
each iteration is very time consuming. Therefore, a fixed 
sensitivity matrix is usually used from the beginning of die 
iteration process, and capacitance A^ is simply calculated by 
multiplying the sensitivity matrix widi die current permittivity 
distribution. Now, equation (29) becomes 

5*^1 =9k- (S^S + yJD-'S^CSA - A). (30) 
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Id literature about ill posed inverse problems, this method 
is called iterative Tikhonov regularization (Hansen 1998). 
Although the Newton-Raphson method is introduced as a non- 
linear method in this section, il should be pointed out that 
iterative Tikhonov Tcgularizationhas actually been used, which 
simplifies image reconstruction as a linear problem. 

3.2, Landmber iteration and steepest descent method 

Landweber iteration is a variation of the steepest gradient 
descent method, which is widely used in optimization dieoiy. 
Suppose that our goal is to minimize ^ |1 5 . ^ - A ||*, i .e. to find 
g which makes the following function minimal. 

/{p) = 5{55-Ar(5p-A) 

= lG7^S^5p -2p^^A + A'^A). (31) 
the gradient of f(g) can be simply calculated as 

'v/to)=S^Sp-S^A = 5^(Sp-A): (32) 

The steepest gradient descent method chooses the direction 
In which f^g) decreases most quickly as the new search 
direction for the next iteration. This direction is opposite to the 
gradient of /(g) at the current point. The iteration procedure 
is therefore 

&+I = & - cti,^f(9k) = P* - atS^iSgk - A) (33) 

where ak is a positive scalar, which decides the ibth step size. 

cr* can be chosen so that V/(&+i) and VfQt) are 
orthogonal. However, choosing an optimal oti needs more 
calculations. A simple way Ls to choose at: as a iixed number 
In ihebegimung of the iteration process. If necessary, its value 
may be changed during iteration (Peng et cd 2000, Liu et ai 
1999). 

Wth a fixed relaxation factor a, the Landweber iterative 
algorithm can be expressed as 



(34) 



The initial guess 50 can be calculated using a simple algorithm, 
e.g. LBP. A problem with Landweber iteration is that its 
convergence property is poor. A simple way to improve its 
conveigence is 10 modify it to give the so-called projected 
Landweber iteration. 



(35) 



where ? is a projection operator. The solution is projected lo 
a convex set after each iteration. 



0 

m 
1 



if/CT)<0 

ifO^/(x)<l 
if/(;c)> L 



(36) 



characteristics. The image error decreases very rapidly in 
the beginning, but after reaching a minimal point, the image 
error increases as the iteration process continues. The optimal 
number of iterations can be determined if there is some prior 
information about the capacitance measurements, but usually 
a fixed number of iterations is chosen empirically, up to a few 
hundred Landweber iteration is the most widely used iterative 
method for ECT (Yang et al 1999, Liu et al 1999). and also 
in masi cases produces the best images as demonstrated in 
section 4. 

The conjugate gradient method, which makes use of 
conjugate direction, is similar to the steepest dc.'icent method. 
It is widely used for solving equations and for optimization. 
Although the conjugate method converges more quickly than 
the stcqpest descent method, it does not give better results 
because it is less regularizing than the steepest descent method. 

5.3, Algebraic reconstruction technique and simultaneous 
iterative reconstruction technique 

The algebraic reconstmction technique (ART) and simultane- 
ous iterative reconstruction technique (SIRT) are commonly 
used for image reconstruction in x-ray computerized tomog- 
raphy (Kak and Slaney 1988). They are simple and effective, 
especially when the system matrix is very large. ART and 
SIKT have been applied to ECT (Reinecke and Mewes 1996, 
Sug/fl/2000). 

ART uses only one set of projecdon data in each iteration 
step. The ART iterative formula for ECT image reconstruction 
based on equation (7) is 



(«*&-! - A*) T 
^ Y *k 



T 



(37) 



This means that the constraint is imposed so that the 
reconstructed image is non-negative and has an upper boundary 
of the material being imaged. Landweber iteration is faster, 
because it only uses the first order derivative. However, 
this docs not necessarily mean that the Landweber method 
converges quickly. Usually, it needs many iteration steps 
before reaching a minimum point. Another drawback 
of the Landweber method rests with its semi-convergence 



where Sk is the Jttfa row vector of the sensitivity matrix. 

The k± iteration of the ART algorithm is based on 
the Jtth row in equation (7). In the k± iteration, only the 
kih normalized capacitance is used to update the normalized 
permittivity distribution vector. If the last row of equation (7) 
is used in the kih iteration, the (^ + 1 )th iteration wiQ be based 
on the first row of equation (7). 

The problem with ART is that it suffers from noise in 
the measurement data. Because only one set of measurement 
data is used in each iteration step, it may not converge if the 
measurement conuins a significant error. SIRT overcomes 
this disadvantage of ART. Assuming that equation (7) is 
composed of M rows, M vectors can be obtained according to 
equation (37) using ART. gk is updated using the average of 
those M vectors. The iterative formula of die SIRT algorithm 

wliere ^ is a relaxation factor, diag(5^) is a vector composed 
of diagonal components of SSF and the division means that 
each numerator is divided by the corresponding denominator. 

If necessary, fi can be replaced by a vector. This means 
tiiat the M vectors based on ART are not treated equally when 
trying to obtain the average vector: Fbr ECT as a non-linear 
system, diis is important for improvement of image quality. 

Actually, equation (38) can be deduced from equation (33) 
by replacing at with a weighting vector. Basically. SIRT is also 
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Figures. Dlnribuiioo models. 

a descent gradient metbod. ItK principle is the same as that of 
Landweber iteration. 

A common issue with the above itetative algorithms is 
that they can only guarantee converging to a local minimum. 
For example, the steepest decent method will give a least 
squares solution, which is unlikely to present a good image, 
If the iteration process is not stopped after certain number of 
iterations. In principle, a suitable solution can be found during 
the iteradve process. The problem is when and how to decide 
the iteradve process should be stopped. This is die key issue for 
most of the iteradve methods for solving ill posed problems. 

3,4. Model-based Uerathn 

A model-based iterative method was proposed by Isaksen 
and Nordtvedt (1993). A few typical models are categorized 
and described by a few parameters. The iteration process is 
carried out by searching for die optimal parameters according 
to the error between the measured capacitance vector and the 
capacitance vector calculated from the cuirem permittivity 
distribution. If some prior information on permittivity 
distribution is known, explicit noodeis can be used. As an 
example, figure 3 shows annular and stradiied distribudons 
in a pipeline. The nK>del for the annular distribudon can be 
described by a single parameter r, and that for the stratified 
distribution by two parameters 9 and d. 

When prior informndon on distribudon is not available, a 
threshold is instead searched for to minimize the discrepancy 
between the measured capacitance and the calculated 
capacitance. This is called implicit model-based iteradve 
image reconstrucdon (Isaksen and Nordtvedt 1993). 

4. Evaluation of algorithms by simulation and 
experiment 

Five image reconstruction algorithms, namely LBP, SVD, 
Tikhonov. iterative Tikhonov and prxyected Landweber 
iteration, were evaluated by both simulation and experiment 
Three evaluation criteria were used: 

(1) relative image error, 

(2) relative capacitance residual and 

(3) correlation coefficient between the test object and the 
reconstruction (Xie et al J 994). 

Because image and capacitance are treated as vectors, their 
norms were used Co calculate die error and die residoal 
respectively. 

Ha - n\\ 

(39) 
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Figure 4. Images reconstructed from simulated data. 

« . 1IA-S*all 
Capacitance residual = ■= — ^^^^ ^" 

Coirelation coefficient = Si^li^' " ~ 



(40) 



(41) 

where g is the true permittivity distribution of the test object, 
^ is the reconsUucted permittivity distribution and g and § ate 
the mean values of ^ and y respectively. 

4,L Evaluation by simulation 

To compare the algoridims, four typical permittivity 
distributions were chosen for simulation with an eight- 
electrode square ECT sensor. Sets of capacitance values for 
all single^electrode combinations were calculated with each of 
these peimittivity distributions, using a finite element method 
(FEM). The lower and higher permittivity materials chosen 
were air and Perspex widi perraitdvity values of 1.0 and 
2.6 respectively. The algorithms were implemented using 
Matiab™ on a PC with a Pentium™ III 750 MHz CPU and 
128 Mbytes memory. 

Figure 4 shows the simulation results. In order to 
compare the speed of the algorithms, die elapsed time was also 
recorded. The image errors, capacitance residuals, correlation 
coefficients and elapsed time of the five algoridiras are listed in 
tables 2-5. For die iterative methods, the number of iterations 
and die relaxation coefKcients are also given. 

Figure 4 shows dtai, as expected, die iterative algoriduns, 
especially die projected Landweber mediod, produce superior 
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Tbble 2. Image cmar (%). 




Single object 


Stratified 


Annular 


TWo objects 


LBP 


87.8 


35.6 


38.6 


85.3 


SVD 


42.4 


1133 


66.8 


57.7 


Hkhonov 


69.2 


46.2 


62.9 


59.8 




0.001 


11 = 0.005 


M = 0.001 


0.0001 


Iterative Tikhonov 


48.5 


20.1 


26.0 


47.8 


(10 iterations, y = 0.01) 










Pjrojected Landweber 


27.2 


14.7 


25.8 


35.6 


{a = 2) 


500 iterations 


400 iieradons 


200 iterations 


500 iterations 



Table 3. Cap&citance residual (%). 





Single object 


Stratified 


Annular 


TWo objects 


LBP 


57.3 


32.1 


29.5 


5J.2 


SVD 










Tikhonov 


19.fi 


28.2 


12 


4.3 




Ai = 0.001 


yx = 0.005 


ti = 0.001 


a = 0.0001 


Iterative Tikhonov 


27.3 


3K0 


19.5 


20.2 


(JO iterations, y =0.0i) 










Projected Landweber 


1.6 


29.7 


19.2 


3.5 


(a = 2) 


500 iterations 


400 iterations 


200 iteratioas 


500 iterations 



'lUile 4. Correlation coefRdent. 





Single object 


Stratiiied 


Annular 


Two objects 


LBP 


0J42 


0.873 


0.861 


0.502 


SVD 


0.897 


0J58 


0.692 


0.789 


Tikhonov 


0.724 


0.790 


0.556 


0.771 




0.001 


M = 0.005 


H = 0.001 


fi = 0.001 


Iterative Tikhonov 


0,912 


0.958 


0.941 


0.882 


(iO iterations, y = 0.01) 










Projected Landweber 


0.960 


0.978 


0.936 


0.925 


(a =2) 


500 iterations 


4(X) iterations 


200 iterations 


500 iterations 


l^le 5. Elapsed time (in seconds). 




Single object 


Stratified 


Annular 


1\vo objects 


LAP 


0.06 


0.06 


0.05 


0.06 


SVD 


0.32 


0J4 


0.33 


0.35 


Tikhonov 


0.33 


0J3 


0.36 


0.32 




ti « 0.001 


fi =5 0.005 


M = 0.001 


/i = 0.001 


Iterative Hkhonov 


26.69 


2636 


26.75 


26.37 


(10 iterations, y = 0.01) 










Projected Landweber 


12.36 


10J8 


5.22 


12.42 


(a =2) 


500 iterations 


400 iterations 


200 iterations 


500 iterations 



interesting to see from table 5 that, although the number of 
iterations for the projected Landweber algorithm is much larger 
than that for the iterative Hkhonov method, it takes a shorter 
lime than the iterative Tikhonov algorithm to converge. 

4,2. Experimental evaluation 

Static experiments were carried out using an ECT system with 
an eight-eiectrode square .sensor ofdimensions 68 cm x 68 cm. 
The signal-to-noise ratio (SNR) of the ECT system is about 
35 dB. Four test distributions were created using Perspex rods 
and Perspex heads: 

( ) ) a Perspex rod of 25 cm in diameter, positioned in the centre 
of the sensor. 



image.s. This is coofinsed by table 2, where in all cases the 
projected Landweber iteration gives the smallest image error. 

Table 3 shows diat SVD produces capacitiince residuals 
very close to zero, although the corresponding image errors as 
listed in table 2 are very significant This indicates that a small 
capacitance error does not necessarily mean a good image. The 
least squares solution, although acceptable mathematically for 
under-determined and ill posed problems, can be physically 
meaningless. Except for SVD, the projected Landweber 
algorithm gives the smallest capacitance errors iu most cases. 

From table 4. it can be seen that in all cases the projected 
Landweber iteration gives the greatest correlation coefficient, 
again indicating the best images. 

Obviously, the elapsed time for iterative algorithms is 
much greater than that for non-iterative algorithms. It Ls 
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Figure S. Images reoonscruacd irom experimental data. 

(2) a Perspex rod of 25 cm in diameter, positioned near the 

sensor wall, 

(3) two Perspex rods of 25 and 15 cm in diameter respectively, 
positioned near the sensor wall, and 

(4) stratified Perspex beads. 

Figure 5 shows the true distributions and the images 
reconstructed from the expoimental data. The capacitance 
eiron were evaluated and arc given in table 6. 

Both the simulation and experimental results show that 
SVD improves the image resolution in the central area, but 
introduces false images elsewhere. Image reconstruction 
based on Tlkhonov regularization strongly depends on the 
regularization param^er. With an appropriate regularization 
parameter, it gives better results than LBP. In general, 
the projected Landweber iteration gives slightly better 
reconstruction than the iterative Hkhonov method. Although 
it requires more iteration steps, it takes less computing time. 

As discussed in section 4.1, a smaD capacitance residual 
does not necessarily imply good image quality. One of the 
reasons is that the capacitance corresponding to the true image 
was simply calculated using a sensitivity matrix multiplied by 
the reconstructed perraitth^ity distribution. This is a typical 
characteristic of ill posed inverse problems. 

4,3. l^ea of errors in capacitance measurements 

The influence of errors in capacitance measuremenK on image 
reconstruction has been investigated, by adding a random 
error into the capacitance vector. A normalized capacitance 



AS • 



.10- 



.OS 



aoo 



' liU0C8pacS8iio&voetDf oonBipondiriptolostphanioni 
' capadtiincBwocBiroBPteniimtodby notsB 
• notsevaclw 




A 



I 



o 

I \ 



10 



25 



Fignre Capadtanoc vecton and noise. 



1 



(tt) 



(4) 



Figure?. Images reconstructed from capacitance with noise. 

(a) LBP; (b) SVD; (c) Tddionov, fi ^ 0,001 ; (d) iterative Tlkhonov 

and (e) projected Landweber. 

vector Xc was calculated using FEM, using the test permiidvity 
distribution shown in figure 4(a). A capacitance vector 
contaminated by noise A, was generated by adding a random 
noise to A^. The relative norm error between and A,, 
which is defined as ||Ae - A,a/||A,||. was 8.9%. Figure 6 
shows the capacitance and noise vectors. Hgure 7 shows the 
images reconstructed using different methods. Tobb 7 gives 
the evaluation criteria for these images. 

Figure 7 and table 7 indicate that image reconstruction 
for BCT is not very srasidve to errors in capacitance 
measurements. The main reason is that the condition number 
of the sensitivity matrix derived from the 32 x 32 grid mesh 
for the eight-electrode ECT sensor used in the simulations is 
only about 50. In this case, the non-linearity dominates image 
reconstruction. Tht error in image reconstruction is caused 
mainly by the linear approximation of sensidvlty manix. 

5. Discussion and future development 

Besides the algorithms categorized above, some other 
techniques have been reponed, e.g. threshold filtering by Xie 
et al (1992) and maximum entropy filtering by Mwambela 
et al (1997). The threshold method was successfully 
used to improve the results based on LBP. In ±c fumre, 
filtering methods including photo-processing methods may be 
combined with the current algorithms to improve die image 
quality (Ge and Yang 2001). Actually, most regularization 
methods may be considered as spectral filters. For example, 
SVD and Tlkhonov regularization according to equation (23) 
may be expressed as 



p 



•VI. 



(42) 
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l^ble 6. Capacitance residual (%) with experimemal data. 



One rod in centre 


One rod near wall 


Two rods 


Stratiiied 


LBP 61.7 


54.5 


44.6 


33.6 


SVD «fO 


^0 


«s0 


«0 


Tikhonov {fX = 0.001) 22.8 


18.7 


17.2 


10.5 


Iterative Hkhonov 49.2 


35.3 


25.4 


44.2 


(10 iterations, y » 0.01) 






43.6 


Projected Laodweber 50.0 


34.8 


31.1 


(a = 2) 200 iteradons 


200 iterations 


50 Iterations 200 iterations 


Table 7. Criteria corresponding to the reconstructed images displayed in figure 7. 




Capacitance 


Correlation 


Elapsed dme 


I mage enor (%) residual (%) 


coefficicm 


(in seconds) 


LBP 87.7 


56.7 


0.540 


0.06 


SVD 58.8 


8.9 


D.813 


a31 


Tikhonov (M = 0.001) 68.9 


193 


0.723 


0.31 


Iterative Hkhonov 48.6. 


26.3 


0.903 


21.3 


O0(terattons.y sO.Ol) 








Projected Landweber 45.1 


7.6 


0.888 


3.5 


(200 iterations, or = 2) 









The equivalent formula in SVD of the Landwebwr iteration 
given in equation (36) is 

9t = - (1 - «^?)*] • ^ • VI (43) 

where k is the iteration step number. A generalized filtering 
solution can be expressed as 

g = f]wi^'!^'Vi (44) 

where forms a spectral filter. 

In Tikhonov regularization, u;i and in Landweber 

iteration^ lU/ = 1 - (1 - a5f)*. For future development in 
image reconstruction for ECT, other spectral filters may be 
considered. The iterative Tikhonov method is ahnost the same 
as Landweber iteration when an appropriate parameter / is 
chosen. 

As discussed in section 4.2, the condition number of 
the sensitivity matrix for an eight-electrode ECT sensor with 
32 X 32 mesh is only about 50. For a 12-electrode ECTT 
sy.stem and the same nieshing strategy, the condition number 
will be about 2000. This means that a 12-eIectrode ECT 
system is inferior to an eight-electrode system in terms of ill- 
conditioning. In the future, further modelling work should be 
carried out so that an optimal sensor structure can be designed* 
including the number and length of electrodes, the meshing 
strategy and the measurement protocol. 

Sdecting the number of electrodes requires a trade- 
off between the number of independent measurements and 
the meastu^meni uncertainty. Fewer electrodes give fewer 
independent measurements, and make the solution more under- 
determined Too many electrodes would result in too small 
sensing area, i.e. too small capacitance to be measured 
accurately. This will increase the uncertainty of capactmnce 
measurements. 

Selecting the length of electrodes also requires a trade- 
off between signal bandwidth and measurement uncertainty. 



Longer electrodes produce an average signal over a greater 
axial length. Shorter electrodes may result in capacitance too 
small to be measured accurately. Typically, tiie electrodes are 
10 cm in length. 

For industrial applications, problems with scaling of 
electrode size must be considered. In principle, the inter- 
electrode capacitance of an ECTF sensor would be independent 
on the sensor diameter if the fringe effect could be ignored. 
This may be the ca.se if the sensor diameter is small, say 
less than 5 cm. However, If the sensor diameter is larger 
than the usual electrode length (10 cm), the fringe effect wUl 
be significant In this case, either long electrodes should be 
used or driven guards should be employed. At present, most 
ECT investigations involve sensors of less than 4" (i.e. 10 cm) 
in diameter although a huge sen.sor of 1 m in diameter was 
used to measure interface levels in an oil separation tank 
(Isaksene/a/ 1994). 

For research investigations, comparatively slow but 
accurate image reconstruction techniques In^lementcd off line 
should be adequate. However, for monitoring and control 
purposes, fast on-line techniques arc desirable. The main aim 
of fiiture work on algorithms will be to provide both more 
accurate and fa.ster image reconstruction algorithms. Parallel 
data processor arrays, e.g. with DSI^, would enable real-time 
iterative image reconstruction. 

Currently, most algorithms are based on a simplified linear 
mathematical model, Becau.se ECT systems are essentially 
non-linear, future work should include investigation into n{)n- 
hnear techniques for both forward problem modelling and 
image reconstruction. Some interesting work has already 
started, such as image reconstmcrion based on electrical field 
lines (Loser et al 2001), and algorithms based on artificial 
ncunU networks (York and Duggan 1997, Nooralahiyan and 
Hpyle 1997, Tapson 1999, Warsiio and Fan 2001). 

To enable ECTT technology to be used in a real industry 
environment, both hardware and software of ECTT systems have 
to be itnproved. On the hardware side, improved SNR and 
increascddataacquisitionratesare requested. OnthesoAwaro 
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side» improved image reconstruction algorithms, eg, rcal-tinie 
Iterative, feature extraction from tomographic images and user- 
friendly >^ndows interface are esKcntia]. 
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